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Introduction 


The  determination  of  potential  coefficients  from  the  analysis  of  satellite 
data,  gravity  data,  or  a combination  of  both  data  types  is  an  active  area  of  research 
in  geodesy.  These  potential  coefficients  are  useful  in  improved  orbit  prediction 
and  estimation;  in  the  computation  of  long  wave  length  features  of  quantities 
dependent  on  the  gravity  field  such  as  geoid  undulations  and  deflections  of  the 
vertical;  and  in  the  information  such  coefficients  imply  about  the  internal  structures 
of  the  earth. 


In  these  determinations  from  two  different  data  types,  satellite  data  and 
terrestrial  gravity  data,  it  is  very  important  to  be  sure  that  we  are  determining 
the  same  quantity  (the  potential  coefficients)  through  consistent  models.  In  this 
paper  we  examine  the  model  problem  for  the  determination  of  potential  coefficients 
from  terrestrial  data.  A recent  discussion  with  numerical  results  for  the  potential 
coefficients  can  be  found  in  Rapp  (1977).  In  that  report  a problem  was  found  with 
certain  corrections  to  be  applied  to  the  terrestrial  anomaly  that  were  needed  to 
determine  potential  coefficients.  To  define  this  further  we  express  the  fully 
normalized  potential  coefficients  (C^a . §£»  ) terms  of  gravity  anomalies  as 
follows  (Ostach  and  Pellinen,  1970): 


(1) 
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4nG(i-l) 


r 

JJ 
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fCosmA'i  _ 

(Ag  + 0!+...)  -i.  r P-a  (sin<P)  da 
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where  G is  an  average  value  of  gravity,  <p  is  the  geocentric  latitude,  P--  are  the 
fully  normalized  associated  Legendre  functions,  Ag  is  the  mean  surface  free  air 
anomaly  given  in  block  da,  and  Gx  is  the  Molodensky  Gx  term.  The  * indicates 
that  the  coefficients  are  given  with  respect  to  the  reference  field  used  for  defining 
the  Ag  values.  Pellinen  (1962)  has  indicated  that  the  neglect  of  the  Gx  term  in 
(1)  can  cause  errors  in  the  low  degree  coefficients  of  10  to  20%.  Numerical  studies 
described  in  Rapp  (1977)  agree  with  this  estimate  but  the  results  are  based  on  a 
number  of  assumptions  relating  the  Gx  term  to  the  terrain  correction. 

In  practice  , the  computation  of  the  Gx  term  for  mean  anomalies  (such 
as  1°  and  5°)  can  be  difficult  in  areas  of  rap  idly  varying  topography.  Consequently, 
it  is  of  interest  to  develop  a procedure  for  the  estimation  of  the  potential  coefficients 
from  terrestrial  data  that  does  not  require  the  computation  of  the  Gx  term,  yet 
retains  the  rigor  implied  by  that  procedure.  In  this  way  we  can  determine  these 
coefficients  in  a manner  that  such  coefficients  can  be  directly  compared  to  satellite 
derived  values.  In  addition,  our  results  should  indicate  a proper  anomaly  (which 
is  easily  computed)  to  use  in  the  combination  of  satellite  and  gravity  data  for 
potential  coefficient  determinations. 
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Procedures 


The  gravitational  potential  of  this  earth  is  usually  expressed  in  terms 
of  the  potential  coefficients  as: 

(2)  V(r,p,A)  = — 1^1  + V ^ (C^B  cos  mA  + sinmA)  Pia  (sintf  )j 

i = 2 1=0 


where  kM  is  the  geocentric  gravitational  constant  and  a is  the  equatorial  radius 
of  an  ellipsoid  approximating  the  geoid.  to  fact,  a is  arbitrary,  but  traditional 
use  has  been  for  it  to  be  taken  as  an  ellipsoid  equatorial  radius.  Equation  (2)  is 
considered  convergent  on  and  outside  the  Brillouin  sphere  R that  just  encloses 
all  the  mass  of  the  earth.  (In  this  discussion  we  will  ignore  the  mass  of  the 
atmosphere.) 


Adopting  a reference  potential  (usually  that  implied  by  an  equipotential 
ellipsoid),  a disturbing  potential,  T,  can  be  found  from  (2).  Using  the  basic 
boundary  condition  of  gravimetric  geodesy  (Heiskanen  and  Moritz,  1967,  p.  86, 
eq.  (2-147c)),  we  find  for  the  gravity  anomaly: 


(3) 


Ag(r,6,A)=^Y  (4-l)(— ) J (Cis  *cosmA  + *sinmA)  P2a  (sin*) 


JL  = a 


i=0 


Equation  (3)  is  regarded  as  formally  convergent  on  and  outside  the  boundary 
sphere.  [It  should  be  noted  here  that  in  practice  the  summation  in  (3)  is  not  to 
= but  to  some  small  finite  degree.  The  resultant  finite  series  has  often  been 
evaluated  at  the  surface  of  the  earth.]  On  the  bounding  sphere  the  gravity 
anomaly  is: 


(4) 


Ag(R,<0,  A)  = §5  Y (£-1) 
R *=3 


a = 0 


We  now  apply  the  usual  orthogonal  relationships  of  spherical  harmonics  to 
derive  the  potential  coefficients  from  anomaly  values  given  on  the  boundary 
sphere.  We  have: 
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This  equation  can  be  compared  to  (1)  and  the  differences  seen.  Let  (Cia , S£3*)a 
be  the  potential  coefficients  computed  from  (1)  setting  & and  higher  terms  to 
zero.  In  addition  we  let: 


(6) 


Ag  (R,<fl,X)  = Ag  (R,«,X)  + 6 


where  Ag  (R,<0 , X)  are  the  anomalies  assumed  to  refer  to  a spherical  approximation 
of  the  geoid  where  R is  the  radius  of  the  sphere.  If  we  now  substitute  (6)  into  (5) 
and  using  (1)  we  have: 


/H''2  + 
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— ; : 6 1 l R.(sinto)  da 
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where  in  the  coefficient  of  the  integral  we  have  set  (kM/R2)  to  an  average  value 
of  gravity  and  (a ,/R)i  has  been  set  to  1.  These  approximations  are  justified  as 
the  6 term  is  expected  to  be  small.  The  coefficient  of  (Cp  , S^.* ) in  (7)  is 
1.0039  at  degree  2 and  reaches  1.0623  at  degree  60  letting  a = 6378140  m, 

R = 6371000  m,  and  R = 6384405  m.  Computations  related  to  (7)  will  be 
described  in  a later  section. 


We  now  turn  to  the  determination  of  the  6 term  introduced  in  equation 
(61.  This  term  may  be  computed  in  several  ways.  An  upward  continuation 
procedure  using  5°  mean  anomalies  has  been  described  by  Rapp  (1969).  However, 
a more  convenient  procedure  is  possibly  the  use  of  least  squares  collocation. 

For  its  application  we  consider  that  we  are  given  a set  of  surface  anomalies 
which  are  to  be  extended  to  the  boundary  sphere.  We  can  represent  this 
computation  in  the  following  general  form: 


(3) 


Ag*  = Cs.s  (Ci*  + Du)"1  [ Ag s] 


where  Ag*  is  the  mean  anomaly  on  the  bounding  sphere,  Ags,  is  the  mean  free-air 
surface  anomaly,  Ct,s  represents  the  covariance  between  Ags  and  Ags  , C<j  is 
the  covariance  between  the  known  surface  anomalies  Ag3  and  D;.  is  the  standard 
deviation  of  the  known  anomalies.  A similar  procedure  has  been  used  in  Rapp 
(1977)  to  estimate  5°  mean  free-air  surface  anomalies  from  lc  x lc  mean  free-air 
anomalies.  However,  in  the  application  of  (8)  we  need  to  use  a spatial  covariance 
function  as  opposed  to  a surface  covariance  function  used  in  Rapp  (1977). 
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In  the  past  few  years  a number  of  spatial  covariance  functions  have 
been  described.  We  can  write: 


(9) 


c (P,Q)  = V cz  P^C0S  ^ 

2 = 3 


where:  C is  the  covariance  between  two  anomalies  located  at  a distance  t?  and 
rq  from  the  center  of  the  earth,  and  separated  by  a spherical  distance  n- . c^  are 
anomaly  degree  variances  and  R3  is  the  radius  of  a sphere  totally  imbedded  within 
the  earth.  This  sphere  is  often  called  the  Bjerhammar  sphere.  Anomaly  covari- 
ance functions  may  also  be  derived  from  disturbing  potential  degree  variances. 

The  infinite  series  in  (9)  may  be  expressed  in  closed  form  provided  a model  is 
chosen  for  the  anomaly  degree  variances.  Examples  of  these  models  and  the 
resultant  closed  expression  for  C (P,Q)  are  described  in  Tscherning  and  Rapp 
(1974).  In  this  report  a subroutine  COVA  was  given  that  computed  covariances 
for  the  following  c^  model: 


(10) 


= A (2  - 1) 
(2-  2)(£-B) 


A = 425.28  mgal3,  B = 24  with  a Rb  value  of  6369780  m.  This  particular 
ts  quite  well  information  that  exists  concerning  point  anomaly  variances, 

. anomaly  variances,  and  anomaly  degree  variances  of  low  degree  as  implied 
irom  satellite  derived  potential  coefficients.  However,  the  computation  of  the 
covariances  needed  in  (8)  requires  an  excessive  amount  of  computer  time  so  that 
we  have  investigated  the  use  of  a simplified  c^  model  that  yields  simplified 
covariance  expressions.  For  such  a case  we  used  the  following  model  suggested 
by  Tscherning: 


(ID 


c2  = 


A'  (1-1) 
2 


which  implies  the  following  anomaly  covariance  function: 


(12)  C (P,Q)  = A's2  (-±  - 1 -2*  J-) 

where 


s3  = R|/(rp  r,); 
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L = [1  - 2st  + s2]~ 


N = 1 + L - st 


t = COS  p 


For  the  application  of  this  model  appropriate  values  of  A'  and  s or  Rs  need  to 
be  chosen.  How  this  was  done  and  results  from  (8)  for  test  cases  are  described 
in  the  next  section. 


Computations  of  the  5 term 


We  first  propose  to  evaluate  (8)  for  5°  mean  anomalies.  In  doing  this 
we  will  consider  three  different  covariance  functions.  The  first  spatial  covariance 
function  is  defined  through  equation  (10)  and  associated  constants  and  is  obtained 
from  subroutine  COVA  given  in  Tscheming  and  Rapp  (1974).  We  then  considered 
the  covariance  function  described  in  (12)  with  values  of  A'  determined  for  two 
different  R3  values.  To  do  this  we  first  defined  Ra  and  evaluated  (12)  with 
p = 0°  and  r?  = r^  = R + 13.405  km.  The  anomaly  variance  from  COVA  was  also 
evaluated  for  this  case  and  the  value  of  A'  computed  so  that  the  anomaly  variance 
at  the  given  r»  = r?  values  would  be  the  same  from  COVA  and  the  simplified 
model  of  equation  (12).  The  two  Rs  values  were  based  on  a value  given  by 
Lauritzen  (1973,  p.  84)  and  on  a value  for  Rs  given  by  Tscheming  (1972,  p.  58). 
We  summarize  in  Table  1 various  constants  to  be  used  in  computations  to  follow. 


Table  1 

Parameters  Associated  with  Covariance  Computations 


Quantity 

Value 

1 1 
Description 

R 

6371000  m 

Mean  Earth  Radius 

R 

6384405  m 

Radius  of  Bounding  Sphere 

Ra 

6369780  m 

Implied  by  Constants  of  COVA 

R 3 

6335960  m 1 

Ra  from  Lauritzen  with  A' 

A' 

16. 95  mgal2  J 

computed  for  use  in  (12) 

Ra 

6348827  m \ 

Ra  from  Tscheming  with  A' 

a' 

10.16  mgal2  / 

computed  for  use  in  (12) 
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For  future  referencing  we  designate  the  covariances  from  (12)  using 

A'  = 16.95  mgar  and  Ra  = 6335960  m as  Model  B,  and  with  A'  = 10. 16  mgar 

and  Rb=  6348827  m as  Model  C. 

For  these  initial  computations  we  chose  4 5°  anomalies  of  upward 
continuation  to  the  bounding  sphere  using  (8).  In  each  case  the  anomaly  on  the 
bounding  sphere  was  computed  using  the  5 closest  known  5°  anomalies  (including 
the  surface  value  itself).  The  mean  anomaly  covariances  needed  in  (8)  were  computed 
by  the  numerical  integration  of  the  point  covariance  function  as  described  in 
Heiskanen  and  Moritz  (1967,  p.  277)  except  our  application  is  to  the  computation 
of  the  covariance  between  mean  anomalies.  This  integration  is  the  most  time 
consuming  aspect  of  the  computation.  Various  tests  were  made  with  different 
grid  spacings  in  this  numerica1  integration  procedure.  For  a given  grid  spacing 
computations  for  a single  5°  block  took  4 times  longer  using  subroutine  COVA 
than  when  the  simplified  model  represented  by  equation  (12).  However,  we  also 
found  that  a smaller  grid  interval  was  required  with  the  use  of  COVA  to  achieve 
a given  accuracy  than  was  required  when  using  the  simplified  model.  Considering 
this  fact,  the  use  of  COVA  in  the  upward  continuation  process  would  take  29  times 
longer  than  the  simplified  model  when  accuracies  on  the  order  of  =0.2  mgals  are 
sought  in  the  process.  An  initial  description  of  the  4 test  blocks  is  given  in  Table  2. 


Table  2 

Description  of  5°  Test  Blocks 


Block 
1 No.* 

| 

Mean  Surface 
FA  Anomaly 

Accuracy 

<P° 

A0 

E i 

| 208 

-18.1  mgals 

-3.4 

45-50 

96-103 

2011  m ; 

310 

- 4.3 

= 4.3 

35-40 

32-  88 

2967 

! 408 

21.5 

= 3.7 

30-35 

319-325 

-3192 

562 

1 

-32. 1 

±3.7 

15-20 

83-  89 

-15S1 

* Rapp  (1977) 


As  can  be  seen  from  Table  2,  two  of  the  blocks  have  a large  mean  elevation, 
while  two  blocks  are  ocean  blocks  with  their  depths  given  . 

Computations  were  made  to  determine  6 as  defined  by  equation  (6)  using 
(8)  with  three  different  covariance  functions.  In  (8)  the  D.j  values  were  set  to 
zero  as  we  wished  to  find  solely  the  effect  of  the  upward  continuation  and  not  a 
filtered  anomaly  value.  The  r?  value  was  determined  by  adding  to  the  mean  earth 
radius  the  average  elevation  of  the  block  (or  zero  if  we  had  an  oceanic  block!. 

The  r^  value  was  taken  as  the  radius  of  the  bounding  sphere  defined  in  Table  1. 
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1 


The  procedure  for  computing  r?  can  be  regarded  as  a spherical  approximation. 

An  improved  ellipsoidal  approximation  can  be  found  by  adding  the  mean  elevation 
of  the  block  to  the  ellipsoidal  geocentric  distance  iplus  geo  id  undulation).  We 
found,  however,  that  the  ellipsoidal  computations  changed  our  final  results  on 
the  order  of  0.2  mgals.  In  the  results  given  in  Table  3 the  spherical  approximation 
has  been  used. 


Table  3 

Estimates  of  5=  Ag  (R,tb,A)  - Ag  (R,<p,A) 
for  4 5°  equal  area  blocks 


Block  No. 

5 (mgals) 

] 

Covar  ance  Used 

208 

0.5 

Original  COVA 

208 

0.6 

Model  B 

0. 5 

Model  C 

310 

0.4 

Original  COVA 

0.4 

Model  B 

310 

0.4 

Model  C 

-1.1 

Original  CO'  A. 

-1.1 

Model  B 

-0.9 

Model  C 

562 

1.6 

Original  COVA 

562 

1.6 

Model  B 

562 

1.3 

Model  C 

I 

It  is  immediately  clear  that  the  results  computed  from  the  three  different  covariance 
functions  differ  only  slightly.  For  computations  dealing  with  the  complete  1654 
anomaly  field  we  choose  to  use  the  constants  of  Model  B defined  previously. 

■ 


From  Table  3 we  see  that  the  influence  of  the  covariances  on  the  upward 
continued  anomaly  is  small.  We  next  examine  two  additional  quantities  that  must 
be  specified  before  a global  computation  can  be  carried  out.  The  first  is  the 
numerical  integration  subdivision  interval  and  the  second  is  the  number  of  known 
5°  anomalies  to  use  in  the  upward  continuation  procedure.  We  first  look  at  the 
integration  interval.  This  arises  when  a point  covariance  function  is  used  to 
compute  the  covariance  between  two  quantities  that  represent  area  averages. 
Consider  the  two  blocks  Ax  and  As  shown  in  Figure  One  for  which  the  covariance 
C»,  is  designed  for  some  unspecified  type  of  variable. 


^1 
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Let  C(P,Q)  be  a point  covariance  function  between  points  P(inAx)  and  Q(inA2). 
Then  we  can  write  (following  from  equation  (7-82)  in  Heiskanen  and  Moritz,  1967): 


If  we  now  let: 

a = ndy 
b = ndx 
a'  = ndy' 
b'  = ndx' 

we  can  express  (13)  in  the  following  numerical  integration  form: 

r.  r.  r.  r. 

(14)  Catas  = r H ) C (x.  y.Xj-y.') 
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The  accuracy  of  this  computation  will  improve  as  n increases,  but  so  will  the 
computation  time.  Our  interest  is  to  choose  the  smallest  n consistent  with  the 
final  accuracies  desired.  To  do  this  we  computed  the  upward  continued  anomalies 
for  the  four  test  blocks  given  in  Table  2 using  different  n values,  and  the  5 
closest  5°  anomalies,  and  the  ellipsoidal  reference  model  (a  = 8373140m, 
f = 1/298.25).  The  results  are  shown  in  Table  4. 


Table  4 

Effect  of  Grid  Size  on  Upward  Continuation  of  5C  Anomalies 


l i n ! 

Block  N'o. 

1 

“t 

6 

3 

10 

203 

-17.27  mgals 

-17.40 

-17.43 

-17.43 

310 

- 3.33 

- 3.90 

- 3.91 

408 

20.49 

20.  64 

20.67 

20.67 

562 

-31.20 

-31. 33 

-31.35 

We  see  from  this  only  small  changes  in  the  upward  continued  anomalies. 
For  our  purposes  we  select  for  future  use  an  n of  5 which  falls  between  two 
tested  values.  Such  a choice  would  appear  to  yield  accuracies  in  this  part  of  the 
upward  continuation  of  about  =0.1  mgals. 

We  next  examined  the  effect  of  the  number  of  known  blocks  to  be  used 
in  the  upward  continuation.  We  carried  out  computations  with  1,5,  and  9 5° 
blocks  closest  to  the  block  being  upward  continued  to  the  bounding  sphere.  These 
results  are  shown  in  Table  5. 


Table  5 


Effect  of  Number  of  Known  Blocks  on 
Upward  Continuation  of  5°  Anomalies 


; Number  of  Known  Blocks 

Block 

1 

5 

9 

1 

208 

-17.  01  mgals 

-17.40 

-17.32 

1 

310 

- 4S.  12 

- 3.90 

- 3.37  ' 

403  | 

20.52 

20.  64 

20.36 

562 

-31. 13 

-31. 33 

-31.34 
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We  see  that  the  difference  between  using  5 and  9 blocks  is  small.  Because  of 
symmetry,  however,  we  choose  to  use  for  further  computations,  the  9 closest 
5°  anomalies  for  the  global  upward  continuation. 

With  the  above  information  available,  the  upward  continuation  of  the 
1654  surface  5°  anomalies  given  in  Rapp (1977)  was  carried  out  under  the  follow  ing 
spec  ifications: 

Covariances:  Model  B 

Grid  Division : 5 

Number  of  Know'n  Blocks.-  9 

In  addition,  the  geocentric  distance  to  the  surface  anomaly  block  was  computed  by 
adding  the  height  of  the  block  to  the  geocentric  radius  of  an  ellipsoid  whose  equatorial 
radius  was  6378140  m with  a flattening  of  1/298.25. 

The  mean  difference  between  the  surface  anomalies  and  the  bounding 
sphere  anomalies  was  0.0  mgals  while  the  root  mean  square  difference  was 
= 0.9  mgals.  The  maximum  difference  was  5. 1 mgals  which  occured  for  a block 
whose  surface  anomaly  was  67  mgals.  A complete  listing  of  the  surface  anomalies 
and  the  bounding  sphere  anomalies  can  be  obtained  from  the  author. 


Application  of  Upward  Continued  Anomalies  to  Potential  Coefficient  Determinations 


We  now  wish  to  use  the  anomalies  on  the  bounding  sphere  to  derive 
potential  coefficients  as  shown  in  equation  (7).  In  doing  this  it  is  convenient  to 
write  (7)  in  the  follow' ing  form: 


(15) 

where 


f-1  * 

f-w0  ' 

rACjt, 

c 

" W-t 

' + 
A 

!as2s 

Y + 

Jt 

W 

= r / R f 

/I 

f -i 

r^i 

• > 

w AS  Jt 

L V R ) 

\ a 

) -1: 

is2 

a J* 

1 

rr  c r 

cos 

ml, 

V 1 

"A5},JU 

4i"G(4- 

i) 

0 < 
Jo  L 

a 

sin 
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Using  the  potential  coefficients  given  by  Rapp  (1977)  and  the  values  of  R,  R, 
and  a (6378140  m)  given  previously  the  above  correction  terms  were  computed 
and  used  to  obtain  the  corrected  coefficients  given  on  the  left  hand  side  of  (15). 
In  order  to  judge  the  magnitude  of  these  correction  terms,  root  mean  square 
values  by  degree,  of  the  total  correction,  correction  (I),  and  correction  (II)  are 
given  in  Table  6.  In  addition,  the  ratio  ( in  percent)  of  these  corrections  to  the 
magnitude  of  the  coefficients  is  shown  in  Figure  2. 

Table  6 

Root  Mean  Square  Value  of  Potential  Coefficient 


Correction  Terms  By  Degree 
(x  10  u) 

i 

Total  Correction 

Correction  (I) 

— 

Correction 

(II) 

2 

717 

796 

1346 

3 

981 

797 

1682 

4 

375 

357 

635 

5 

470 

234 

605 

6 

234 

266 

443 

7 

312 

200 

472 

8 

151 

113 

214 

9 

222 

147 

344 

10 

114 

103 

195 

11 

132 

93 

206 

12 

93 

53 

132 

13 

91 

66 

140 

14 

59 

50 

99 

15 

82 

57 

132 

16 

71 

54 

117 

17 

66 

47 

104 

18 

66 

42 

100 

19 

52 

39 

SI 

20 

59 

31 

86 

21 

57 

34 

35 

22 

56 

35 

85 

23 

52 

30 

78 

24 

51 

28 

74 

i 

25 

54 

29 

79 

26 

58 

27 

80 

27 

55 

25 

74 

28 

60 

27 

32 
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29 

56 

30 

82 

30 

56 

26 

76 

31 

50 

22 

67 

32 

53 

27 

75 

33 

50 

24 

70 

34 

54 

27 

77 

35 

44 

25 

65 

36 

51 

27 

74 

37 

42 

24 

61 

38 

46 

26 

67 

39 

42 

22 

60 

40 

43 

26 

64 

From  Figure  2 we  see  that  the  correction  terms,  relative  to  the 
coefficients,  increase  with  degree.  The  first  correction  term  shows  almost  a 
linear  dependence  on  l while  correction  term  (II)  shows  a more  rapid  increase 
than  (I).  The  total  correction  is  less  (on  the  average)  than  correction  term  (II) 
as  there  is  some  cancellation  between  the  two  correction  terms.  We  see  then 
that  the  correction  terms  have  larger  effects  at  high  degrees  than  at  the  lower 
degrees.  This  is  opposite  to  what  was  found  in  Rapp  (1977)  when  techniques  for 
estimating  the  Gx  term  in  equation  (1)  were  investigated.  The  results  of  this 
paper  seem  to  be  more  realistic. 

We  should  point  out  here  that  the  correction  terms  are  small,  not  only 
with  respect  to  the  potential  coefficients,  but  also  with  respect  to  the  estimated 
standard  deviations  of  the  potential  coefficients  as  determined  from  the  standard 
deviations  of  the  anomalies  used  in  the  computation.  Such  comparisons,  in  terms 
of  %,  are  given  in  Table  7. 


Table  7 

Total  Correction  Term,  and  Standard  Deviation  of, 
Potential  Coefficients  Expressed  as  a Percent  of  the  Coefficients 


l 

} 

Correction 

Std.  Dev. 

l 

Correction 

Std.  Dev. 

5 

1 

1.8% 

13.7% 

25 

5.  5% 

80.5% 

10 

1.6 

29.1 

30 

3.4 

86.0 

15 

2.7 

45.5 

35 

7.4 

87.2 

20 

4.6 

77.7 

40 

7.5 

83.4 

1 
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It  is  clear  from  Table  7 that  the  correction  terms  are  considerably  smaller  than 
the  uncertainty  in  the  coefficients.  On  the  other  hand,  the  corrections  computed 
here  are  20  times  larger  than  the  corrections  to  the  potential  coefficients  arising 
from  the  use  of  a spherical  reference  surface,  and  the  neglect  of  the  atmosphere 
(Rapp,  1977). 

Our  next  computation  involved  the  comparison  of  the  potential  coefficients 
computed  from  5°  gravity  data  (Rapp,  1977)  and  those  coefficients  corrected  using 
equation  (15),  with  the  GEM  9 potential  coefficients  (Lerch  et  al. , 1977).  The 
GEM  9 coefficients  are  complete  to  degree  20  and  are  derived  solely  on  the  basis 
of  satellite  observations.  These  comparisons,  in  terms  of  root  mean  square 
potential  coefficient  differences,  are  given  in  Table  8. 

Table  8 

Root  Mean  Square  Potential  Coefficient  Differences 
Between  the  GEM  9 Coefficients  and  Those  Coefficients 
Implied  by  5°  Mean  Anomalies  Given  on  the  Earth's 
Surface  and  On  the  Bounding  Sphere 


1 

Surface  Anomalies 
Used 

Bounding  Sphere 
Anomalies  Used 

2 

. 2934  x 10“s 

.2946  x 10"s 

3 

.2625 

.2638 

4 

. 1612 

. 1619 

5 

. 1473 

. 1494 

6 

.0935 

.0926 

7 

.0574 

.0573 

8 

.0592 

.0593 

9 

.0612 

.0604 

10 

.0444 

.0444 

11 

.0332 

.0331 

12 

.0254 

.0252 

13 

.0222 

.0222 

14 

.0252 

.0251 

15 

.0228 

.0224 

16 

.0205 

.0202 

17 

.0217 

.0214 

18 

.0172 

.0169 

19 

.0145 

.0143 

20 

.0164 

.0161 
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From  this  table  we  see  a slight  degradation  in  the  agreement  when  the  correction 
terms  are  applied  at  the  lower  degrees.  However,  from  degree  6 and  beyond,  the 
corrected  coefficients  agree  with  the  GEM  9 coefficients  better  than  the  uncorrected 
coefficients.  However,  the  improvement  is  small  and  could  be  considered  negligable 
when  the  accuracy  of  the  potential  coefficient  determinations  is  considered.  Never- 
theless the  results  are  encouraging. 


Application  of  Upward  Continued  Anomalies  to  Anomaly  Comparisons  with 
Anomalies  Implied  by  Potential  Coefficients 


For  a number  of  years  potential  coefficient  solutions  have  been  judged  by 
computing  anomalies  from  such  coefficients  and  comparing  such  anomalies  with 
surface  anomalies.  In  this  computation,  equation  (3)  or  a spherical  approximation 
has  been  used  to  compute  the  anomalies  from  the  coefficients.  Such  a procedure 
has  been  done  ignoring  the  convergence  question  as  it  has  been  considered  to  have 
negligable  influence  on  the  final  results.  We  now  turn  our  attention  to  such 
comparisons  using  the  results  of  this  paper. 

To  do  this  we  first  computed  anomalies  in  5°  blocks  as  implied  by  the 
GEM  9 coefficients  to  degree  20  using  equation  (3)  where  r was  evaluated  as  the 
geocentric  distance  to  the  block  whose  anomaly  was  being  determined.  Let  this 
anomaly  be  designated  Ags*  . Next  compute  the  5°  anomaly  on  the  bounding  sphere 
using  equation  (4).  Using  this  value,  compute  from  equation  (6)  the  anomaly  on 
the  surface  using  the  6 value  implied  by  the  upward  continuation  process.  This 
anomaly,  which  we  will  designate  Agsa  , can  be  rigerously  compared  to  the 
surface  anomaly  values,  as  given  in  Rapp  (1977).  Such  comparisons  between  the 
observed  surface  values  Ags  , Ags*  and  Agse  are  shown  in  Table  9 in  terms  of 
root  mean  square  differences  (when  the  mean  difference  has  been  excluded). 


Table  9 

Root  Mean  Square  Difference  Between  Potential 
Coefficient  Anomalies  and  Surface  5°  Anomalies 


1 

m(max)  (mgals) 

I 

5 

10 

15 

E ((  Ags  - Ags*)s) 

107 

109 

123  ; 

| E ((  Ags-  Agsa)3) 

88 

91 

101 
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In  this  table  m is  the  maximum  standard  error  of  the  terrestrial  anomaly  used 
in  the  comparison.  The  number  of  such  anomalies  is  1003  (m-5  mgals); 

1343  (ms  10  mgals)  and  1433  (m  £ 15  mgals). 

It  is  clear  that  the  better  agreement  between  the  anomalies  occurs 
when  the  potential  coefficients  are  used  to  evaluate  anomalies  on  the  bounding 
sphere  which  are  then  downward  continued  to  the  surface,  instead  of  just 
evaluating  the  anomalies  on  the  surface  using  an  expression  of  questionable 
converging  properties. 


Summary  and  Conclusions 


This  paper  has  been  written  to  show  how  potential  coefficients  can  be 
determined  from  terrestrial  gravity  data  rigorously  considering  the  topography 
by  upward  continuing  the  surface  free-air  anomalies  to  a sphere  that  completely 
surrounds  the  earth.  This  upward  continuation  can  be  performed  by  least  squares 
collocation.  We  found  that  the  use  of  sophisticated  covariance  functions  was  too 
expensive  in  terms  of  computer  costs  so  that  it  was  necessary  to  choose 
simplified  equations.  However,  tests  showed  that  results  from  the  simplified 
equations  agreed  quite  well  with  values  obtained  from  the  more  realistic 
covariance  functions.  Using  these  simplified  models,  a set  of  1654  5°  equal 
area  anomalies  on  the  bounding  sphere  have  been  computed.  The  maximum 
difference  between  the  surface  anomaly  and  the  anomaly  on  the  bounding  sphere 
was  5.1  mgals  with  the  root  mean  square  difference  being  =0.9  mgals. 

The  correction  term  5,  and  the  potential  coefficients  implied  by  the 
surface  data,  were  used  in  equation  (15)  to  derive  correction  terms  to  the  surface 
data  potential  coefficients.  The  correction  terms  were  small,  being  about  2 ~ 
of  the  actual  coefficients  at  degree  5 rising  in  a roughly  linear  way  (percentage 
wise)  to  about  7%  at  degree  40. 

The  corrected  coefficients  were  compared  to  the  GEM  9 potential 
coefficients,  as  were  the  coefficients  implied  by  the  uncorrected  surface 
anomalies.  We  found  that  above  degree  5,  the  corrected  coefficients  agreed 
slightly  better  with  GEM  9 than  did  the  coefficients  implied  by  the  uncorrected 
surface  anomaly  data. 

A final  test  was  made  by  computing  surface  anomalies  directly  from  a 
spherical  harmonic  expansion  evaluated  on  the  surface  of  the  earth  (and  thus  of 
questionable  accuracy),  and  also  by  first  computing  the  anomalies  on  the  bounding 
sphere  and  subsequently  downward  continuing  them  to  the  earth's  surface  using 
the  5 term  derived  through  collocation.  These  two  anomaly  sets  were  then 
compared  to  the  surface  anomalies  given  in  Rapp  (19  77).  We  found  that  the  better 
agreement  was  obtained  when  the  corrected  bounding  sphere  anomalies  were  used 
(mean  square  difference  =109  mgalsi,  than  when  directly  computed  anomalies  on  th 
surface  from  potential  coefficients  were  used  (mean  square  difference  = 91  mgals i 


1 


This  paper  has  shown  that  the  problems  raised  by  the  convergence 
question  are  real  at  the  = 1 mgal  level  for  5°  anomalies.  Thus  for  the  most 
accurate  determination  of  potential  coefficients  from  surface  gravity  material, 
equation  (7)  or  (15)  should  be  used. 

Finally,  an  important  implication  of  these  studies  exists  for  those 
computations  that  combine  surface  gravity  and  satellite  data  for  potential 
coefficient  determinations.  In  most  current  applications,  equation  (2)  is  used 
to  relate  anomalies  and  potential  coefficients.  However,  this  approach  is  in 
error  with  the  more  correct  procedure  being  the  use  of  equation  (4)  with  the 
needed  anomaly  being  found  from  equation  (6). 
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